Rabi Oscillations in Systems with Small Anharmonicity 
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When a two-level quantum system is irradiated with a microwave signal, in resonance with the 
energy difference between the levels, it starts Rabi oscillation between those states. If there are 
other states close, in energy, to the first two, the Rabi signal will also induce transition to those. 
Here, we study the probability of transition to the third state, in a three-level system, while a Rabi 
oscillation between the first two states is performed. We investigate the effect of pulse shaping on 
the probability and suggest methods for optimizing pulse shapes to reduce transition probability. 
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I. INTRODUCTION 

Most qubits (i.e. basic elements in a quantum com- 
puter) are not true two-level systems. Yet, only the first 
two energy states are commonly considered relevant for 
quantum computation. As a result, any transition to the 
upper levels during the gate operations is a leakage of in- 
formation outside the computational space, and therefore 
a source of error. 

One of the common methods to perform gate opera- 
tions in a qubit is via Rabi oscillations 0. The speed 
of operation is determined by the Rabi frequency £Ir, 
which is proportional to the amplitude of the applied 
microwave signal. Rabi oscillations have been observed 
in many quantum systems, including superconducting 
qubits 0, y, 0, IE El, excitons in single quantum dots 
0, @ , and very recently single electron spins in nitrogen- 
vacancy defect centers in diamond 0- 

In a multi-level quantum system, Rabi oscillations may 
not be limited to only the first two states. For example, in 
a harmonic oscillator, with equally spaced energy eigen- 
values, applying a Rabi signal in resonance with the level 
spacings will occupy many states. When the system is 
strongly anharmonic, on the other hand, i.e. when the 
third state is far above the first two, the probability of 
transition will be vanishingly small. 

To have a quantitative measure of anharmonicity, we 
define an anharmonicity coefficient by 



(E 2 i — E w )/E 



in. 



(1) 



where, 25y — Ei — Ej, with Eq being the ground state 
and E^q, the i-th excited state energy. 5 is zero for a 
harmonic oscillator and — > + oo for an ideal two level 
system. 

Not every qubit realization has large S. For example, 
in a current biased Josephson junction qubit Q|, E21 is 
always smaller than E10 leading to a negative S close 
to zero. Charge-phase qubits also suffer from small an- 
harmonicity, merely because of operating in the charge- 
phase regime; for the "quantronium" qubit of Vion et 
al. 5 « 0.2, and for the flux based charge-phase qubit 
of Ref. 0, a (5 = 0(1) was suggested. 

The purpose of this paper is to study how much small- 
ness of 5 can affect transition to the upper state, and 
how it can be prevented. We study the problem in a 
three-state quantum system with small anharmonicity. 



In Sec. II, we perform analytical calculations using Ro- 
tating Wave Approximation (RWA). Section III, goes be- 
yond RWA using numerical methods. The effect of pulse 
shaping on the transition probabilities is addressed in 
Sec. IV. Section V, discusses practical examples within 
superconducting qubit implementations. A brief sum- 
mary together with some concluding remarks are pro- 
vided in Sec. VI. 



II. ANALYTICAL CALCULATION 

Let us consider a quantum system with three states 
i — a,b,c, irradiated with a microwave signal in 
resonance with the energy difference between the first 
two levels. The Hamiltonian of the system is written as 

[E c >E b >E a =0) 



H = Eb\b)(b\ + E c \c)(c\ + V(t) 



(2) 



where Vlt) = V c 



oe 



-iujot _ 



c.c. is the microwave signal 



(K = 1). Writing the wave function as — c a (t)\a) 
Cb{t)\b) + c c (t)\c), the equations of motion for Cj are 



ic a = V ab c b , 

ic b = V* b c a + E b c b + V bc c c 
ic c = V bc c b + E c c c , 



(3) 



where Vy(t) = (i\V(t)\j). We have taken V ac = 0; the 
transition probability will be small anyway because of 
large frequency difference. For simplicity, we write E b = 
u>o and E c = (2 + 5)u>q. In this section, we assume 6 <C 1 
to ensure small anharmonicity. 

i2w t an( j wr jt e 



Let us define c b 
ue~ 



c b e^\ c c 



V a b/oJo = ue~ iuJot + c.c. and V bc /uj = ve~ iu ° t + c.c. Us- 
ing RWA, i.e. ignoring the fast oscillating terms, we find 

d T c a = -iuc bl 

d T c b = —iu*c a — ivc c , 

d T c c — —iv*c b — iSc c , (4) 

where r = u>ot. The equation of motion for c b can be 
extracted from 

[dl + ml + (M 2 + \v\ 2 )d T + is\ u \ 2 ] ~ Cb = 0. (5) 



Writing c b = ke lXT ', x needs to satisfy 

x 3 ~ Sx 2 - (\u\ 2 + \v\ 2 )x + 5\u\ 2 



(6) 



2 



General solutions are 
x n = — 1 5 + 2z cos 
where 



(2n-l)J]}, n = 1,2,3 (7) 



z = v / 3(|"| 2 + |f| 2 )+<5 2 , 



1 (96(\u\ 2 -\v\ 2 /2) -S 3 

= - arccos — - 

3 \ z A 



To find the coefficients, let us write 

3 



(8) 
(9) 



Cb 



EL p -ix n r 

n=l 
1 3 
u ^— ' 

n— 1 

7/ 7 1 * * 



(10) 



n=l 



which satisfy @- Assuming that the system starts from 
the ground state, we impose the initial conditions: c a = 1 
and q, = c c = 0, which yield 

3 3 3 

E fc„ = 0, ^i„fc„ = 0, E a^fcn = |u| 2 (11) 

n—1 n—1 n—1 

Solving these equations for k n , we find 

M 2 + X 2 X 3 



ki 



(Xl - X 2 )(xi - £3) 



(12) 



k 2 and k 3 can be obtained using the permutation 
l->2->3->l. 

Let us write c c = ^a n e~ lXnT , where a n = (a: 2 — 
\u\ 2 )k n /uv. The probability of finding the system in the 
upper state is 



Pc(r) = |c c | 2 



where 



\n=l 



< Prr 



(13) 



(14) 



determines an upper bound for P c (r). We first study the 
solution is some special cases. 



A. Case I, 5 = 

This is the simplest case that the problem can be 
solved. From (0-©, we find 



xi = 0, x 2 ,3 = T^fl, 



(15) 



where D,r = \/\u\ 2 + \v\ 2 . These can also be found easily 
from © directly. Using (|12|1 . we find ki = \v/£Ir\ 2 and 
k 2 = k 3 = (l/2)\u/Q R \ 2 . As a result 



Ml 



(H 2 + \u\ 2 cosn R T), 



-i— sinS! fi T, 



n 2 



■(1 — cos £Irt). 



(16) 



The system oscillates with only one frequency £Ir. The 
probability of finding the system in the upper state |c c | 2 
can become large: P max = 4|uv| 2 /£l^,. This is expected 
in a system with zero anharmonicity. 



B. Case II, v = 

Using together with 

cos36> = 4cos 3 6»-3cos6>, (17) 

we find cos# = —5/z, which immediately gives 

xi = — \u\, x 2 = 6 : x 3 = \u\. (18) 

These could also be found directly from For k's, we 
get: fci = k 3 = 1/2 and ^2 = 0, leading to 

c a = cosQrt, 



u 

~Cc = 0. 



(19) 



The results show usual Rabi oscillation between the first 
two states with frequency Q r = \u\. The probability 
of finding the system in the upper state is always zero 
(P c = 0), as expected because v = 0. 



C. Case III, S > u, v 

In the regime u, v <C 6 -C 1, one can find asymptotic 
solutions. A systematic expansion in u/S and v/S gives 



xi = \u\ 1 



\v\ 2 \ \v\ 2 



2S 2 J 26 



x 2 = -\u\ 1 



2S 2 J 26 



x * = 5 ( 1 + ^ 



Leading to the Rabi frequency 



n R = \u\ i - 



26 2 



(20) 



(21) 
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FIG. 1: P max vs 5/u for different values of v/u. The curves 
are symmetric with respect to 5 — > —5. 

The dependence of the Rabi frequency on the amplitude 
of the microwave signal now has the form 

Sl R oc Vo (1 - /3V 2 ) , (22) 

where the coefficient /3 depends on the details of the sys- 
tem. The deviation from the proportionality relation is 
a signature of transition to the upper states. Such a de- 
viation has been experimentally observed recently in a 
current biased dc-SQUID structure [l2|. 

The probability of finding the system in the upper state 
is given by 

P c « M- sin 2 n RT . (23) 

It oscillates with the Rabi frequency The maximum 
probability 

-Pmax «-^« 7 ^_j (24) 

occurs at half a Rabi period r = 7t/Qh, where is the 
largest. This is not the case for small 5 (see e.g. case I). 
Here, 7 = |w/w| 2 is a constant depending on the details 
of the Hamiltonian. In most physical systems \v\ ~ \u\ 
and therefore 7 = 0(1). 

D. General case 

It is not straightforward to find a closed analytical so- 
lution for the general case. Instead we plot the results 
for Pmax, calculated using Q-© together with l|12fl and 
1)14)1. Figure ^ shows P max as a function of 5/u with dif- 
ferent values of v/u. At small v/u, the curves are peaked 
near 5 = u, while for larger v/u the peak appears near 
(5 = 0. In all cases P ma x becomes very small at large 5/u, 
as expected. 

III. NUMERICAL CALCULATION 

In this section we calculate the quantum evolution of 
the system numerically using density matrix approach. 




FIG. 2: Probability P a (solid), P b = p 22 (dashed), and P c = 
P33 (dot-dashed), as a function of time. The parameters are 
u = v = 0.1, 5 = 0.1 (a) and S = 0.5 (b). For clarity, P c in 
(b) is magnified by a factor of 10. 

This allows us to study the system beyond RWA and/or 
at large 5. The dynamics of the 3x3 density matrix p is 
described by 

if t =[H,p]. (25) 
We integrate this equation starting from 

/l 0\ 
po = , (26) 
\ / 

which describes the system at the lowest energy state. 
Probabilities of finding the system in different states are 
given by: P a = p n , Pb = P22, and P c = p 33 . Figure El 
displays the time evolution of these probabilities. The 
fast oscillations are the effect of high frequency terms, 
which were ignored in the previous section due to RWA. 
Figure^ shows the Rabi oscillation when 5 = 0.1. After 
(almost) half a Rabi period, significant amount of the 
probability goes to the third state. By increasing 5 to 
0.5, the probability of finding the system in the upper 
state is significantly reduced (Fig. Eb; the curve in the 
figure is magnified for clarity). 

The maximum probability of the system in the upper 
state is given by P max = Max r [P c ]. Figure 01 shows the 
dependence of P max on 5. The solid lines are analyti- 
cal curves using (|14p. and the dashed ones represent the 
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results of numerical calculations. While the two curves 
coincide at small 5, they soon deviate from each other as 
S increases. However, the overall behavior of the curves, 
especially the asymptotic P max ~ H 2 /5 2 dependence re- 
mains unchanged even at large S. To emphasize on this 
aspect, we have plotted P ma x<5 2 /|^| 2 vs S in Fig. 0] for 
different values of parameters. All the curves overlap at 
large 5 suggesting P max ~ \v\ 2 /S 2 ~ (n R /6) 2 , in agree- 
ment with (|24|l : the coefficient 7, however, is now a slow 
function of the parameters, but still 0(1). 




0.5 1 1.5 2 
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FIG. 3: Pmax vs 8 for different values of u and v. Solid 
(dashed) curves are analytical (numerical) results. 



u=0.1, v=0.1 

u=0.1, v=0.02 



u=0.02, v=0.02 

u=0.02, v=0.1 




0.5 1 1.5 2 



5 

FIG. 4: P m ax5 2 /M 2 vs 8 for different values of u and v. 



IV. EFFECT OF PULSE SHAPE 

So far we have assumed that the microwave signal 
starts at r = 0, and continues forever. To perform a gate 
operation, however, one needs to apply the Rabi signal 
for only a short duration of time. In that respect, our 
calculation can only describe hard pulses, in which the 
microwave switches on and off abruptly. The probability 
P c then oscillates with the pulse duration at the Rabi 
frequency. The maximum probability usually happens in 
the case of a 7r-rotation, i.e when the probability is maxi- 
mally transferred to \b). A hard pulse, however, is neither 
practical, nor the best pulse shape, as was indicated in 



Ref.lH Indeed, by using other types of pulses, the prob- 
ability of transition to the upper level, at the end of the 
process, can be significantly reduced. Among a few pulse 
shapes examined in [l^, Gaussian pulses demonstrated 
the most promise. To understand the role of pulse shap- 
ing, let us compare the effect of a Gaussian pulse on the 
probability P c , with that of a hard pulse, for the case of 
a 7r-rotation. 

To enforce a Gaussian envelope for the microwave sig- 
nal V(t), we write 

u(t) = f {aT/r w )e-^-^) 2 ^ for < r < r p 
\ otherwise 

where t p and t w are the duration and width of the pulse 
respectively, T is the total angle of rotation in the Bloch 
sphere (e.g. V = tt for a 7r-rotation) , and a is a normal- 
ization constant. 




100 200 300 400 500 

X 



FIG. 5: Probability of the third level as a function of time for 
a hard and a Gaussian 7r-pulse. 

Figure shows the probability P c as a function of time 
for a Gaussian and a hard pulse, both of which having 
the same duration and resulting in a 7r-rotation (r = tt) 
at the end of the pulse. In our numerical calculation 
we take v = u, 6 = 0.05, r p — 500, t w = r p /6, and 
a = 0.398. These numbers correspond to the optimal 
pulse shape suggested in |13| . The maximum of P c for 
the Gaussian pulse, happens slightly after the center of 
the pulse, while in the case of the hard pulse, it occurs 
near the end. Although the maximum is larger for the 
Gaussian pulse, the probability P c f at the end of the pro- 
cess is much smaller. Orders of magnitude reduction of 
the final probability can be achieved using such a tech- 
nique. 

In Ref. Il3l t w was fixed (to t p /6 or t p /4) and r p was 
varied to minimize P c f. A r p f=a 87r/|5| was shown to 
provide the first minimum with shortest duration. Alter- 
natively, one can fix r p and find a t w which gives mini- 
mum P c f- This may work better for shorter pulses. For 
example, for r p = 100, 5 = 0.1, and v = u, a Gaussian 
pulse with t w = Tp/6 gives P c f — 0.093, while the mini- 
mum probability P c f — 0.0026 is achieved at r w — 0.31r p 
and a = 0.467. Such a pulse shape starts and ends with 



5 



jumps (see Fig. |HJ), but still gives smaller P c at the end 
of the process. 
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FIG. 6: Pulse shapes optimized for a 7r-rotation with 8 = 
0.1, t p = 100. 

A Gaussian pulse shape is not the optimal pulse shape 
for minimizing P c f. One can design other pulses with 
more free parameters to achieve a smaller probability. 
To have some idea about how small can P c f be made by 
appropriately shaping the pulse, we defined an arbitrary 
pulse by the scries 



u(t) = (r/Tp) 



A„ cos(27rnr / V p ) 



(27) 



Keeping only the first two terms in the series, (using 
the same conditions as above: t p = 100, 6 = 0.1, and 
v = u) one can already reach a probability as small as 
P cf = 1.2 x 10~ 5 with Ai = -0.3833 and A 2 = 0.1293 (see 
Fig. EJ) . With N = 33 terms in the series, the probability 
was reduced to 2.4 x 10~ 6 . The resulting pulse shape, 
shown in Fig. [7| is complicated and may not be useful 
experimentally. It should also be emphasized that with 
the pulse shape of l|27|l . there is not a unique minimum 
for P c f. Depending on the starting point and the method 
of minimization, one may fall into a local minimum with 
complicated pulse shape. 




FIG. 7: Pulse shape of Eq. ffij . optimized with N = 33. 

Here, we only considered the case of V — ir. For 
quantum operations, other pulses may also be required. 



It is not just enough to change the amplitude of the 
pulse, keeping its shape and duration, to obtain opti- 
mized pulses with other T's. Indeed, for each type of op- 
eration, one needs to design a specific pulse shape that 
provides minimum P c f. 



V. DISCUSSION 

In a practical quantum computer, the maximum num- 
ber of operations is limited by the decoherence time of the 
qubits as well as the speed of operations. It is generally 
believed that if ~ 10 4 operations can be performed within 
the decoherence time, quantum computation can con- 
tinue indefinitely with the help of quantum error correc- 
tion algorithms. A parameter that is commonly quoted 
as a measure for the maximum number of operations is 
the quality factor of the qubits, usually defined as 



1 



(28) 



where t v is the dephasing time of the qubit (in units of 
1/uo). Q v , however, is related to only one type of single 
qubit operations, namely phase rotation. Other neces- 
sary operations such as single qubit state flip or multi- 
qubit gate operations are usually much slower. Even for 
the phase rotation, the extent to which one can control 
the rotation, i.e. change -Eio, may be much smaller than 
the rotation frequency itself. 

The singlequbit state flip can be performed u sing Rabi 
oscillations U 0,111 or non-adiabatic evolution |ll| The 
latter is fast (~ u>o), but requires large anharmonicity 
to avoid unwanted Landau-Zener transition to the upper 
states. Rabi oscillations, on the other hand, are much 
slower, but can be used in small anharmonicity systems. 
It is possible to define a quality factor for the Rabi oscil- 
lations the same way as Q v was defined in l(2"%|) 



Qr — -^p' RTR ~ ®-b.Q<pi 



(29) 



where tr is the Rabi decay time which is typically the 
same order as t v . 

In an ideal two level system, Qr is limited by the 
maximum allowed amplitude of the microwave signal (re- 
stricted by RWA and/or experimental limitations). Usu- 
ally an Qr as large as 0.1 or even larger is conceiv- 
able. In practical systems, especially those with small 
anharmonicity, however, increasing the microwave power 
will cause transition to the upper states as we discussed. 
Therefore Qr is limited by how much probability of the 
upper levels can be tolerated. If we restrict P max to 
- 10~ 4 , then gives fl R ~ 10~ 2 <S. Therefore to 

achieve fJ_R ~ 0.1 (Qr ~ 0.1Q V ), we need a 6 > 10. 
Such a large anharmonicity cannot be supported by many 
qubit implementations (see below for a few examples). 

Using a shaped (instead of hard) pulse can significantly 
reduce the final P c . To define a quality factor similar to 
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(|29|l . we use the fact that in the case of a hard pulse, a ir- 
rotation is implemented when = tt/t p . We therefore 
define 

Qshapcd = - f — ) t v = ( — ) Q v . (30) 
2 \ t pJ \ t pJ 

Therefore a Qshapod = 0.1Q y requires a pulse with dura- 
tion t p — 10n ss 30 for a 7r-rotation. It was shown in |13| . 
that a Gaussian pulse with t w = t p /6 provides minimum 
P c with shortest time if r p « 87r/|<5] . A quality factor of 
O.IQ^ is therefore achievable in a system with 8 « 0.8. 
Other pulse shapes may provide better performance at 
smaller 5, as was discussed before. Below, we provide a 
few examples among superconducting qubits. 

In the current biased Josephson junction qubit of 
Ref. |j, the energy differences are wio ~ 6.9 GHz and 
^21 ~ 6.28 GHz^ leading to 8 w —0.09. Also, one can 
easily justify jl^l that \v\ = V2\u\ ~ \u\, as expected. 
For a hard pulse, requiring P max ~ 10~ 4 and using l(2^|) 
(with 7 w 1), one finds Q,r ~ 10~ 3 wo, which is extremely 
slow. The quality factor Qr will also be very small 
(~ 10~ 3 <3 V ). Aiming for a larger quality factor, one can 
make use of shaped pulses. A Gaussian pulse with du- 
ration t p — 100 (Q s hapcd ~ 0.03(5^) and with optimized 
width (t w = 36.4) gives P c / = 0.0073, which may not 
be small enough. The pulse shape of Eq. l|77|l . optimized 
with only first two components (Ai = —0.2331, A2 = 
0.2916), on the other hand, gives a probability as small 
as P c f = 1.6 x 10~ 5 , for the same pulse duration. It is 
not easy to reach a small P c f with a shorter pulse. 

In the charge-phase (quantronium) qubit of Ref. 0, 
S 0.2, fl R ~ 100 MHz, and w » 16 GHz. We there- 
fore obtain \u\ w 0.0063, and with |«| ~ \u\, using i)24|) we 
find P max ~ 5 x 10~ 4 for a hard pulse, which is reason- 
ably small. The quality factor for the Rabi oscillation, 
however, is Qr w 150 much smaller than Q v = 25000 
quoted in Q . Increasing the Rabi frequency will increase 
the probability P max - With the help of a Gaussian pulse 
shape (with optimal width t w — 15.3), a pulse duration 
of Tp = 50 (quality factor Qshapcd ~ 0.06Q V ) is achiev- 
able with P c f = 0.0026. Again, significant improvement 
in the probability (P c / = 9.3 x 10~ 6 ) can be achieved 
using Eq. I|27|) . optimized keeping only two components 
in the series (Ai = -0.4058, A 2 = 0.1241). 

In practice, the shape of the pulse should be motivated 
experimentally. For example, the jumps at the ends of 
the pulses shown in Fig. [S] can only be realized approxi- 
mately. Such limitations should be considered as a con- 
straint in the optimization process. The minimization 



procedure may also be preformed experimentally; trying 
different pulses with a few free parameters and probing 
the transition probability to the upper levels. 



VI. SUMMARY AND CONCLUSIONS 

We have performed analytical and numerical investi- 
gations of Rabi oscillations in a three level system. We 
showed that the probability P c of finding the system in 
the upper level oscillates with the Rabi frequency f^. 
The maximum probability P max happens close to half a 
Rabi period. We demonstrated that P max ~ (Qr/5) 2 , 
even beyond RWA and when S is large. 

We also studied the effect of pulse shaping on P c . We 
showed that with an appropriate pulse shape, one can 
achieve small probability P c at the end of the process, 
although in the middle of the operation it may become 
large. The duration and shape of the pulse can be opti- 
mized to obtain smallest P c / in a shortest time. For each 
type of necessary operation, a specific pulse shape should 
be designed. In any case, smallness of 5 limits how short 
the pulse can be and therefore affects the speed of qubit 
operations. 

It is also necessary to take into account the effect of 
decoherence on the studied phenomenon. In practice, 
however, only a few Rabi oscillations happen during the 
operation. Thus, as long as the decoherence time of the 
system is much longer than the Rabi period, our conclu- 
sions remain valid even in the presence of decoherence. 

In this article, we only considered three levels. If the 
anharmonicity of the system is very small, one needs to 
consider more than three states. In Ref. [liq . ~ 10 states 
were taken into account in the numerical simulations. 
Finally, we should mention that having a multi-level, in- 
stead of two-level, quantum system is not necessarily a 
disadvantage, as long as coherent control of all the levels 
is possible. There have been proposals to use multi-level 
systems for quantum computation |l4j . 
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